Biostar:课程13、14
翻译:生物女博士
注:本系列课程翻译已获得原作者授权。
#为作者的注释
#*为译者的注释
Lecture 13 - 安装和使用bwa
#从SourceForge上下载bwa
# 从SourceForge获取的URL是比较复杂的,因为SourceForge希望你去访问网址并从网址上下载,这样你可以看到网页投放的广告。
# 你可以访问网址,并手动下载。
#*最新版本是bwa-0.7.15
#*响应时间可能会比较久
cd ~/src
curl -OL http://downloads.sourceforge.net/project/bio-bwa/bwa-0.7.15.tar.bz2
# 这个安装包是用bzip2压缩的,所以需要加-j 来解压。
tar jxvf bwa.0.7.15.tar.bz2
# 我们需要编译一下程序
cd bwa-0.7.15/
# 用make命令来编译。
make
# 把bwa程序链接到bin目录下。
ln -s ~/src/bwa-0.7.15/bwa ~/bin
# 运行一下程序,确保其可用。
mkdir -p ~/edu/lec13; cd ~/edu/lec13
# 运行这个工具,看看它是怎样用的。
bwa
# 手册不错呦!
man ~/src/bwa-0.7.15/bwa.1
#*点击"q"可以退出。
# 创建一个总目录来存储参考序列。
mkdir -p ~/refs/852
# 获取1999年大爆发的埃博拉病毒的参考基因组序列。
#*貌似最近用efetch工具下载NCBI上的数据有问题。大家可以手动下载,然后放到目的文件夹。
efetch -db nucleotide -id NC_002549 -format fasta > ~/refs/852/ebola-1999.fa
# 对基因组创建索引,这样bwa就可以进行基因组搜索。
bwa index ~/refs/852/ebola-1999.fa
# 来看一下运行结果。
ls ~/refs/852/
# 使之也能通过blast进行搜索。
makeblastdb -in ~/refs/852/ebola-1999.fa -dbtype nucl -parse_seqids
# 现在,来看看发生了什么。
ls ~/refs/852/
# 从一个文件中衍生了16个文件,我们最好还是把它们放在一个独立文件夹。
ls -l ~/refs/852/ | wc -l
# 获取埃博拉基因组序列的第一行。
head -2 ~/refs/852/ebola-1999.fa > query.fa
# 用bwa-mem命令将其map回它的基因组。
bwa mem ~/refs/852/ebola-1999.fa query.fa > results.sam
# 与该序列的blast过程比较一下
blastn -db ~/refs/852/ebola-1999.fa -query query.fa > results.blastn
Lecture 14 - 了解SAM格式
# 从lecture13结束的地方继续讲。
head -2 ~/refs/852/ebola-1999.fa > query.fa
# 做序列比对。
bwa mem ~/refs/852/ebola-1999.fa query.fa > results.sam
# 注意,bwa还是会在屏幕上输出“流”(stream)。
#* stream:这里是指软件运行过程中屏幕上不断出现的一些文字。
# 有两种输出流。它们被称为:标准输出(standard out)和标准错误输出(standard error)。
# 下边的代码是告诉你,你可以如何重定向这两个输出到文件中。
bwa mem ~/refs/852/ebola-1999.fa query.fa 1> results.sam 2> errors.txt
# 有很多技巧可以将一个流重定向到文件,尽管在绝大多数情况下,这是不必要的。但有时确实需要存储错误消息。
# http://www.tldp.org/LDP/abs/html/io-redirection.html
# 根据你有的序列创建“查询序列”。
# 添加或删除一些碱基来创造“插入”和“缺失”。
# 做序列比对。
bwa mem ~/refs/852/ebola-1999.fa query.fa > results.sam
# 查看某些列。
# 了解每列的含义很重要。
#*这里作者讲的过于简略,推荐大家看参考阅读1和2或者官网。
# (这三列是)query id, 比对标识和目标id
#*query id其实就是reads序列的名字。
cat results.sam | cut -f 1,2,3
# 比对起始位置,mapping质量, 简要比对信息表达式
cat results.sam | cut -f 4,5,6
# 双端reads信息:名字,位置和模板长度
cat results.sam | cut -f 7,8,9
# Query序列, query质量, 其他属性
#*Query序列(也就是reads序列)
cat results.sam | cut -f 10,11,12,13,14
# 如果你列两个文件给bwa-mem,它会以配对模式来比对它们。
# 从课程网站上获取示例数据集。
curl -O http://www.personal.psu.edu/iua1/courses/files/2014/data-14.tar.gz
# 解压和比对这个数据集,根据结果回答作业上的问题
bwa mem ~/refs/852/ebola-1999.fa read1.fq read2.fq > results.sam
参考阅读1:http://www.cnblogs.com/emanlee/p/5366610.html
参考阅读2:http://blog.sina.com.cn/s/blog_670445240101l30k.html
欢迎关注我们